{smcl}
{com}{sf}{ul off}{txt}{.-}
      name:  {res}<unnamed>
       {txt}log:  {res}R:\WSV2\TBu_AKe\Spatial_NEW\spatial_3_wmat_beam.smcl
  {txt}log type:  {res}smcl
 {txt}opened on:  {res}29 Mar 2023, 12:23:44
{txt}
{com}. 
. 
. cd "R:\WSV2\TBu_AKe\Spatial_NEW"
{res}R:\WSV2\TBu_AKe\Spatial_NEW
{txt}
{com}. 
. 
. use "C:\Users\hy65byfe\Desktop\smerge_0712\spatial_2_collapse_s.dta" // data set with 4950 cells, created in spatial_2_*  
{txt}
{com}. 
. 
. scalar p = 0 // if 0, straight distance without decay : x^0 = 1, so all are equal.
{txt}
{com}. scalar sigma = 47*0.59 // sigma at comp. line 
{txt}
{com}. scalar kappa = 51.7*0.59 // kappa at comp. line 
{txt}
{com}. scalar D = 3 // width of increment 
{txt}
{com}. 
. ***************************
. **** reduce dimensions ****
. ***************************
. 
. keep if year == 2011 & ccode == 1 // reduce data in memory 
{txt}(549,450 observations deleted)

{com}.   
. sort cell // now sorted: first bunching (1-50), then restricted (51-1015) then rest. 
{txt}
{com}. 
. gsort -bunching +cell 
{txt}
{com}. gen bcell = _n if bunching == 1 // identify B
{txt}(4,900 missing values generated)

{com}. 
. gsort -restricted +cell 
{txt}
{com}. gen rcell = _n if restricted == 1 // identify R
{txt}(2,100 missing values generated)

{com}. 
. 
. **********************************************
. *** define upper and lower bounds of cells ***
. **********************************************
.                 
.         ** gen mid-point for E 
. gen a_mid = a_bin 
{txt}
{com}. gen e_mid = `=kappa' + `=sigma'*a_mid //  
{txt}
{com}. 
.         ** gen upper and lower level of E at mid-point (e.g. 6kg, 6.5kg ...)
. replace e_bin = k_bin + a_bin*`=sigma' // gravity center: upper edge of mid-point 
{txt}(4,950 real changes made)

{com}. gen e_ml = k_bin + a_bin*`=sigma' - `=D' // lower edge of mid-point 
{txt}
{com}. 
. 
. replace eei_bin = (e_bin/(51.7+47*a_bin))* 100 // assert that bins are consistent with eei boundaries. 
{txt}(4,944 real changes made)

{com}. 
.         ** label as safety net 
. label variable e_bin "E_i at midpoint"
{txt}
{com}. 
. label variable a_min "a_i left"
{txt}
{com}. label variable a_max "a_i right"
{txt}
{com}. label variable k_bin "kappa (std. E)"
{txt}
{com}. label variable e_mid "E_m"
{txt}
{com}. 
. 
. keep eei_bin a_* e_* k_bin s_bin eei_bin cell bcell rcell // cut 
{txt}
{com}. 
. sort cell // order is crucial 
{txt}
{com}. 
. ***********************************
. ******* start loops ***************
. ***********************************
. 
. **********************************       
. *  3a: neighborhood slope / fan *
. **********************************
. 
. ***********************************
. *** gradients from unit circle  ***
. ***********************************
.         
.         ** sw 
. scalar g0 = tan((225*(_pi/180)))
{txt}
{com}. scalar g1 = tan((202.5*(_pi/180)))
{txt}
{com}.         ** horizontal 
. scalar g2 = tan((180*(_pi/180)))
{txt}
{com}.         ** mixed (nw) 
. scalar g3 = tan((157.5*(_pi/180)))
{txt}
{com}. scalar g4 = tan((135*(_pi/180)))
{txt}
{com}. scalar g5 = tan((112.5*(_pi/180)))
{txt}
{com}.         ** vertical 
. scalar g6 = -tan((90*(_pi/180))) // has to be minus due to quadrant switch at 90. 'approach' from left on unit circle.
{txt}
{com}. scalar g7 = tan((67.5*(_pi/180)))
{txt}
{com}. scalar g8 = tan((45*(_pi/180)))
{txt}
{com}. 
. forvalues i = 0/8 {c -(}
{txt}  2{com}.         scalar delta`i'= (1/g`i')
{txt}  3{com}. {c )-}       
{txt}
{com}. 
. * scalar list 
. 
. **************************************
. * 3b: define scalar for every cell ***
. **************************************
. 
. 
. sum cell, meanonly 
{txt}
{com}. scalar cmax=`r(max)'
{txt}
{com}. local c = `r(max)' // scalar for each cell.  
{txt}
{com}.         
. sum bcell, meanonly 
{txt}
{com}. scalar bmax=`r(max)'
{txt}
{com}. 
. sum rcell, meanonly 
{txt}
{com}. scalar rmax=`r(max)'
{txt}
{com}. 
.         *** scale to normalize slope of compliance line to 1. From formula: EEI = (E / k + ba ) * 100
. scalar sd_a = 100/`=sigma'
{txt}
{com}. scalar sd_e = 100
{txt}
{com}. 
. 
. forvalues i = 1/`=cmax' {c -(}
{txt}  2{com}.         * 5-digit cell-id
.         scalar cell`i' = 10000 + `i'  
{txt}  3{com}.         
.         * capacity 
.         quietly sum a_bin if cell == `i'
{txt}  4{com}.         local a = r(mean)
{txt}  5{com}.         scalar a`i' = `a'/sd_a
{txt}  6{com}.         
.         quietly sum a_min if cell == `i'
{txt}  7{com}.         local a = r(mean)
{txt}  8{com}.         scalar a_min`i' = `a'/sd_a
{txt}  9{com}.         
.         quietly sum a_max if cell == `i'
{txt} 10{com}.         local a = r(mean)
{txt} 11{com}.         scalar a_max`i' = `a'/sd_a
{txt} 12{com}.         
.         * std.  
.         quietly sum s_bin if cell == `i'
{txt} 13{com}.         local s = r(mean)
{txt} 14{com}.         scalar s`i' = `s'
{txt} 15{com}.         
.         quietly sum k_bin if cell == `i'
{txt} 16{com}.         local s = r(mean)
{txt} 17{com}.         scalar k`i' = `s'
{txt} 18{com}.         
.         * eei 
.         quietly sum eei_bin if cell == `i'
{txt} 19{com}.         local x = r(mean)
{txt} 20{com}.         scalar x`i' = `x'
{txt} 21{com}.         
.         * e [four courners and midpoint]
.         quietly sum e_bin if cell == `i'
{txt} 22{com}.         local e = r(mean)
{txt} 23{com}.         scalar e`i' = `e'/sd_e
{txt} 24{com}.         
.         quietly sum e_low_left if cell == `i'
{txt} 25{com}.         local e = r(mean)
{txt} 26{com}.         scalar e_ll`i' = `e'/sd_e
{txt} 27{com}.         
.         quietly sum e_up_right if cell == `i'
{txt} 28{com}.         local e = r(mean)
{txt} 29{com}.         scalar e_ur`i' = `e'/sd_e
{txt} 30{com}.         
.         quietly sum e_up_left if cell == `i'
{txt} 31{com}.         local e = r(mean)
{txt} 32{com}.         scalar e_ul`i' = `e'/sd_e
{txt} 33{com}.         
.         quietly sum e_low_right if cell == `i'
{txt} 34{com}.         local e = r(mean)
{txt} 35{com}.         scalar e_lr`i' = `e'/sd_e
{txt} 36{com}.         
.         quietly sum e_mid if cell == `i'
{txt} 37{com}.         local e = r(mean)
{txt} 38{com}.         scalar e_mid`i' = `e'/sd_e
{txt} 39{com}. {c )-}
{txt}
{com}. 
. 
. ***********************************
. **** 3c: define neighborhood ******      
. ***********************************
.                  
. ****************************     
.         * inclusion conditions 
.         
. forvalues k = 1/`=bmax' {c -(}
{txt}  2{com}.      quietly sum cell if bcell==`k'
{txt}  3{com}.      local m = r(mean)
{txt}  4{com}.           * m contains cellid of cell 1  
.           forvalues l = 1/`=cmax' {c -(}
{txt}  5{com}.           quietly sum cell if cell==`l'
{txt}  6{com}.           local n = r(mean)
{txt}  7{com}.                  * n contains cellid of cell 2 
.                   forvalues d = 0(1)8 {c -(}
{txt}  8{com}.                         ** qualifying conditions: two parallel beams from edges of bunching cell touching the compliance line 
.                   scalar qr`d'_`m'_`n' = (e_ur`m'+((a`n'-a`m')/delta`d')) // uppper-right boundary 
{txt}  9{com}.                   scalar ql`d'_`m'_`n' = (e_ul`m'+((a`n'-a`m')/delta`d')) // upperleft boundary 
{txt} 10{com}.           {c )-}
{txt} 11{com}.      {c )-}
{txt} 12{com}. {c )-}        
{txt}
{com}. 
. 
. ******************************* 
. ***  4a: defining neighbors. ***
. *******************************
. 
. 
. forvalues k = 1/`=bmax' {c -(}
{txt}  2{com}.      quietly sum cell if bcell==`k'
{txt}  3{com}.      local i = r(mean)
{txt}  4{com}.           * i contains cellid of cell i in B     
.           forvalues l = 1/`=cmax' {c -(}
{txt}  5{com}.           quietly sum cell if cell==`l'
{txt}  6{com}.           local j = r(mean)
{txt}  7{com}.                   * j contains cellid of cell j in R 
.                   
.                   ** knock-out condition: use only relevant segment (moved up  for efficiency: saves 7*4950 scalars to compute)
.                   scalar x_`i'_`j'              = cond( x`j' > x`i'                             &  inrange(x`j', 59, 87)                ,       1,0) // only neighbors in B. --> hadamard. 
{txt}  8{com}.                   
.                   forvalues d = 0(1)8 {c -(}   
{txt}  9{com}.                    ** test boundaries against conditions : neighbor if gravity center falls in beam   
.                    scalar s`d'_`i'_`j'  = cond(ql`d'_`i'_`j' <= e`j'            &  qr`d'_`i'_`j' >= e`j'                ,       1,0)  // 
{txt} 10{com}.                    
.           {c )-}
{txt} 11{com}.      {c )-}  
{txt} 12{com}. {c )-}
{txt}
{com}. 
. 
. *********************************
. *** 4b: transform to matrix *****
. *********************************
. display cmax // check dimension 
{res}4950
{txt}
{com}. 
. *** neighborhood matrices
. mat define X  = I(`=cmax') // knock out criterium
{txt}
{com}. mat define S1 = I(`=cmax') // sw 1   
{txt}
{com}. mat define S2 = I(`=cmax') // sw 2
{txt}
{com}. mat define S3 = I(`=cmax') // nw 1 
{txt}
{com}. mat define S4 = I(`=cmax') // nw 2
{txt}
{com}. mat define S5 = I(`=cmax') // nw 3
{txt}
{com}. mat define S6 = I(`=cmax') // nw 4 
{txt}
{com}. mat define S7 = I(`=cmax') // no 1
{txt}
{com}. mat define S8 = I(`=cmax') // no 1
{txt}
{com}. 
. ***************************
. **  5: fill matrix ********
. ***************************
. 
. forvalues k = 1/`=bmax' {c -(}
{txt}  2{com}.      quietly sum cell if bcell==`k' & bcell != . 
{txt}  3{com}.      local i = r(mean)
{txt}  4{com}.      forvalues l = 1/`=cmax' {c -(}
{txt}  5{com}.           quietly sum cell if cell==`l' 
{txt}  6{com}.           local j = r(mean)
{txt}  7{com}.                   ** run through combinations and put in distance **
.                   mat X[`k',`l']=(x_`i'_`j')
{txt}  8{com}.                   mat S1[`k',`l']=(s1_`i'_`j')
{txt}  9{com}.                   mat S2[`k',`l']=(s2_`i'_`j')
{txt} 10{com}.                   mat S3[`k',`l']=(s3_`i'_`j')
{txt} 11{com}.                   mat S4[`k',`l']=(s4_`i'_`j')
{txt} 12{com}.                   mat S5[`k',`l']=(s5_`i'_`j')
{txt} 13{com}.                   mat S6[`k',`l']=(s6_`i'_`j')
{txt} 14{com}.                   mat S7[`k',`l']=(s7_`i'_`j')
{txt} 15{com}.                   mat S8[`k',`l']=(s8_`i'_`j')
{txt} 16{com}.          {c )-}
{txt} 17{com}.          forvalues c = 1/`=cmax' {c -(}
{txt} 18{com}.         ** replace own-cross with zero (identical cell)
.         mat X[`c',`c']=0
{txt} 19{com}.         mat S1[`c',`c']=0
{txt} 20{com}.         mat S2[`c',`c']=0
{txt} 21{com}.         mat S3[`c',`c']=0
{txt} 22{com}.         mat S4[`c',`c']=0
{txt} 23{com}.         mat S5[`c',`c']=0
{txt} 24{com}.         mat S6[`c',`c']=0
{txt} 25{com}.         mat S7[`c',`c']=0
{txt} 26{com}.         mat S8[`c',`c']=0
{txt} 27{com}.         {c )-}
{txt} 28{com}.     {c )-} 
{txt}
{com}. 
.         
.         ** HADAMARD for each sector : eliminates those that are not from R or do not meet condition. 
. forvalues i = 1/8 {c -(}
{txt}  2{com}.         matrix NS`i'=hadamard(X, S`i')
{txt}  3{com}. {c )-}       
{txt}
{com}.         
. ***********************************************
. ****** 6: post results to weight matrix. ******
. ***********************************************
. cd "$store"
{res}R:\WSV2\TBu_AKe\Spatial_NEW\Data
{txt}
{com}. 
.         ********************
.         ** Sectors 1-4 *****
.         ********************
. preserve
{txt}
{com}. 
. svmat NS1, names(s1_)
{txt}
{com}. svmat NS2, names(s2_)
{txt}
{com}. svmat NS3, names(s3_)
{txt}
{com}. svmat NS4, names(s4_)
{txt}
{com}. 
. 
. 
. save spatial_3_wmat_beam_1234, replace
{txt}{p 0 4 2}
(file {bf}
spatial_3_wmat_beam_1234.dta{rm}
not found)
{p_end}
{p 0 4 2}
file {bf}
spatial_3_wmat_beam_1234.dta{rm}
saved
{p_end}

{com}. 
. restore
{txt}
{com}. 
. 
.         ********************
.         ** Sectors 5-8, Hadamard ******
.         ********************
. preserve
{txt}
{com}. 
. svmat NS5, names(s5_)
{txt}
{com}. svmat NS6, names(s6_)
{txt}
{com}. svmat NS7, names(s7_)
{txt}
{com}. svmat NS8, names(s8_)
{txt}
{com}. 
.         ** check ** 
. sum s5_82 s6_82 s7_82 s8_82 

{txt}    Variable {c |}        Obs        Mean    Std. dev.       Min        Max
{hline 13}{c +}{hline 57}
{space 7}s5_82 {c |}{res}      4,950     .000404    .0200987          0          1
{txt}{space 7}s6_82 {c |}{res}      4,950     .000404    .0200987          0          1
{txt}{space 7}s7_82 {c |}{res}      4,950     .000404    .0200987          0          1
{txt}{space 7}s8_82 {c |}{res}      4,950     .010101     .100005          0          1
{txt}
{com}. tab cell if s5_82 > 0

       {txt}cell {c |}      Freq.     Percent        Cum.
{hline 12}{c +}{hline 35}
          7 {c |}{res}          1       50.00       50.00
{txt}         32 {c |}{res}          1       50.00      100.00
{txt}{hline 12}{c +}{hline 35}
      Total {c |}{res}          2      100.00
{txt}
{com}. tab cell if s6_82 > 0 

       {txt}cell {c |}      Freq.     Percent        Cum.
{hline 12}{c +}{hline 35}
          7 {c |}{res}          1       50.00       50.00
{txt}         32 {c |}{res}          1       50.00      100.00
{txt}{hline 12}{c +}{hline 35}
      Total {c |}{res}          2      100.00
{txt}
{com}. tab cell if s7_82 > 0 

       {txt}cell {c |}      Freq.     Percent        Cum.
{hline 12}{c +}{hline 35}
          7 {c |}{res}          1       50.00       50.00
{txt}         32 {c |}{res}          1       50.00      100.00
{txt}{hline 12}{c +}{hline 35}
      Total {c |}{res}          2      100.00
{txt}
{com}. tab cell if s8_82 > 0 

       {txt}cell {c |}      Freq.     Percent        Cum.
{hline 12}{c +}{hline 35}
          1 {c |}{res}          1        2.00        2.00
{txt}          2 {c |}{res}          1        2.00        4.00
{txt}          3 {c |}{res}          1        2.00        6.00
{txt}          4 {c |}{res}          1        2.00        8.00
{txt}          5 {c |}{res}          1        2.00       10.00
{txt}          6 {c |}{res}          1        2.00       12.00
{txt}          7 {c |}{res}          1        2.00       14.00
{txt}          8 {c |}{res}          1        2.00       16.00
{txt}          9 {c |}{res}          1        2.00       18.00
{txt}         10 {c |}{res}          1        2.00       20.00
{txt}         11 {c |}{res}          1        2.00       22.00
{txt}         12 {c |}{res}          1        2.00       24.00
{txt}         13 {c |}{res}          1        2.00       26.00
{txt}         14 {c |}{res}          1        2.00       28.00
{txt}         15 {c |}{res}          1        2.00       30.00
{txt}         16 {c |}{res}          1        2.00       32.00
{txt}         17 {c |}{res}          1        2.00       34.00
{txt}         18 {c |}{res}          1        2.00       36.00
{txt}         19 {c |}{res}          1        2.00       38.00
{txt}         20 {c |}{res}          1        2.00       40.00
{txt}         21 {c |}{res}          1        2.00       42.00
{txt}         22 {c |}{res}          1        2.00       44.00
{txt}         23 {c |}{res}          1        2.00       46.00
{txt}         24 {c |}{res}          1        2.00       48.00
{txt}         25 {c |}{res}          1        2.00       50.00
{txt}         26 {c |}{res}          1        2.00       52.00
{txt}         27 {c |}{res}          1        2.00       54.00
{txt}         28 {c |}{res}          1        2.00       56.00
{txt}         29 {c |}{res}          1        2.00       58.00
{txt}         30 {c |}{res}          1        2.00       60.00
{txt}         31 {c |}{res}          1        2.00       62.00
{txt}         32 {c |}{res}          1        2.00       64.00
{txt}         33 {c |}{res}          1        2.00       66.00
{txt}         34 {c |}{res}          1        2.00       68.00
{txt}         35 {c |}{res}          1        2.00       70.00
{txt}         36 {c |}{res}          1        2.00       72.00
{txt}         37 {c |}{res}          1        2.00       74.00
{txt}         38 {c |}{res}          1        2.00       76.00
{txt}         39 {c |}{res}          1        2.00       78.00
{txt}         40 {c |}{res}          1        2.00       80.00
{txt}         41 {c |}{res}          1        2.00       82.00
{txt}         42 {c |}{res}          1        2.00       84.00
{txt}         43 {c |}{res}          1        2.00       86.00
{txt}         44 {c |}{res}          1        2.00       88.00
{txt}         45 {c |}{res}          1        2.00       90.00
{txt}         46 {c |}{res}          1        2.00       92.00
{txt}         47 {c |}{res}          1        2.00       94.00
{txt}         48 {c |}{res}          1        2.00       96.00
{txt}         49 {c |}{res}          1        2.00       98.00
{txt}         50 {c |}{res}          1        2.00      100.00
{txt}{hline 12}{c +}{hline 35}
      Total {c |}{res}         50      100.00
{txt}
{com}. 
. save spatial_3_wmat_beam_5678, replace
{txt}{p 0 4 2}
(file {bf}
spatial_3_wmat_beam_5678.dta{rm}
not found)
{p_end}
{p 0 4 2}
file {bf}
spatial_3_wmat_beam_5678.dta{rm}
saved
{p_end}

{com}. 
. restore
{txt}
{com}. 
. log close 
      {txt}name:  {res}<unnamed>
       {txt}log:  {res}R:\WSV2\TBu_AKe\Spatial_NEW\spatial_3_wmat_beam.smcl
  {txt}log type:  {res}smcl
 {txt}closed on:  {res}29 Mar 2023, 12:54:51
{txt}{.-}
{smcl}
{txt}{sf}{ul off}